8 Zooarchaeology
8.1 Case studies
The following map shows the sites under investigation, divided by chronology. Please select the desired chronology (or chronologies) from the legend on the right. Legend: R = Roman, LR = Late Roman, EMA = Early Middle Ages, Ma = 11th c. onwards
The faunal dataset used in this study is both extensive and diverse, containing over 466 records. While NISP is a useful proxy for historical animal farming, it is not without its limitations, as previously discussed in the methods section. What is crucial to emphasize here, however, is the presence of overdispersion within the data, which necessitates more sophisticated and nuanced approaches than simply calculating overall means for each animal species. Overdispersion is a common occurrence when analyzing datasets of this nature, given the unique factors at play in each context, including historical and depositional influences. Nevertheless, data modeling requires simplification and causal reasoning, and the best approach is to start with straightforward models that can account for overdispersion and generate credible distributions. To this end, Bayesian hierarchical models were developed for each chronology, context type, macroregion, and geography. As a further step, an additional analysis was conducted that solely focused on altitude and chronology as predictors. By examining the specific contributions of altitude and chronology in predicting animal farming/consumption patterns, this analysis provides valuable insights that complement the earlier models. These findings underscore the importance of considering multiple factors when studying the probability of occurrence of farmed and wild animals in historical contexts. They also demonstrate the potential benefits of simplified models to focus on key predictors. Moreover, several attempts were made to create a coherent understanding of the likelihood of economically valuable animals occurring during the first millennium, in order to provide a more comprehensive perspective on the role of animal farming in shaping historical societies.
8.2 Data exploration
As previously mentioned, the faunal dataset used in this study exhibits overdispersion, which is a common issue when analysing datasets of this type. The term ‘overdispersion’ refers to the presence of greater variability in the data than expected based on a normal curve (Figure 8.1). In this section, I will present the distribution of the animals of interest to provide a visual representation of the dispersion within the dataset. By examining the distribution curves, we can better understand the variability that exist within the dataset, and can use this information to develop more accurate models. In Figure 8.2, it is evident that the distribution of animal remains in the faunal dataset is not symmetrical and does not conform to a normal curve. This non-normal distribution indicates that the standard measures of central tendency, such as mean and median, may not accurately capture the overall distribution of the data. Traditional statistical measures, such as measures of central tendency and dispersion, are commonly used in frequentist approaches to analyse data. However, these measures may not be appropriate for analysing the complex patterns of animal farming and consumption in the faunal dataset due to the presence of overdispersion. As an alternative, Bayesian multilevel models can account for overdispersion by incorporating appropriate probability distributions, such as the betabinomial distribution.
8.3 Chronology
8.3.1 Trends by century
The proposed model uses a betabinomial distribution to model overdispersion in the data and estimates the precision (or shape) parameter \(\phi\) in the Beta distribution. The \(A\) on the left side of the formula is the outcome variable, which represents the animal NISP count for each observed sample \(i\). The model is intercept-only, meaning that there are no further predictors included in the model. The average probability of success \(\bar{p}_{i}\) is modelled through a logit link function, and is assumed to equal the intercept \(\alpha\). The model is also indexed with \({[Century]}\) so that separate estimates can be obtained for each century (1st BCE to 11th CE).
\[ A_{i} \sim BetaBinomial(NISP_{i}, \bar{p}_{i} , \phi_{i}) \]
\[ logit(\bar{p}_{i}) = \alpha_{[Century]} \]
\[ \alpha_{[Century]} \sim Normal(0,1.5) \]
\[ \phi_{[Century]} \sim Exponential(1) + 2 \] The prior distributions selected for the model are weakly informative. Specifically, the prior chosen for the intercept \(\alpha\) is a normal distribution with a mean of 0 and a standard deviation of 1.5. Given that \(logit(\bar{p}) = \alpha\) is modelled using the logit link, a mean of 0 corresponds to a probability of 0.50 on the logit scale. With a standard deviation of 1.5, the prior is flat and only slightly informative.
As for the shape parameter \(\phi\), the prior was selected to ensure that the parameter was at least 2, which corresponds to a flat distribution. This value (2) was then added to an exponential distribution with rate 1 to produce a weakly informative prior.
A prior predictive simulation (Figure 8.3) shows how when no information is provided to the model, counts on the extremes of the binomial histogram are more likely, while all other counts have similar probabilities. Bimodal priors are common practice with Beta distributions as they are non-informative enough. However, updating a bimodal prior even with one observation rapidly changes the posterior.
The model for pigs’ NISP counts shows a positive trend above the millennium average (0.36) for predicted high density intervals (HDIs) from the 1st century BCE to the 3rd century CE. The probability intervals exhibit a steady decrease that begins in the 3rd century and continues until the 8th century, when the mean of the probability distribution coincides with the first millennium average. Despite positive trends reported in the 8th and 9th centuries, interpretation of these values requires caution since the credible intervals are wider due to fewer observed samples dated to these centuries. The HDIs mean for the 10th and 11th centuries is again around the millennium average, indicating a decrease in the probability of pig occurrence.
Cattle, on the other hand, exhibit different trends, with a HDI interval in the 1st century BCE above the millennium average (0.20), although with a large credible interval, and a negative trend from the 1st to the 3rd century CE. HDIs are again increasingly positive between the 4th and 7th centuries, but decrease in the 8th century. As for pigs, HDIs means are close to the millennium average in the 10th and 11th centuries.
The modelled probability of occurrence of caprine reveals a positive trend in the 1st century BCE, followed by a negative trend (below the average of 0.28) from the 1st to the 3rd century CE. The probability of occurrence then increases around the 4th/5th centuries up to the 6th/7th centuries, with HDIs stable around the mean from the 8th to the 11th century.
Edible wild mammals and domestic fowl predicted probabilities are relatively stable trends throughout the analysed period, with HDI means consistently below 0.10. Domestic fowl HDIs remain below the millennium mean of 0.06 from the 1st century BCE to the 2nd century CE, and then increase steadily from the 3rd century to remain around the mean until the 7th century. From the 7th to the 11th century, the HDIs remain positive, but with a much larger credible interval. Similarly, wild mammals’ HDIs are positive in the 1st century BCE and remain stable around the millennium average (0.05) until the 6th century CE, after which they also show a slight increase until the 11th century.
In addition to plotting the probabilities of occurrence extracted by the models, the distributions of the precision parameter \(\phi\) are plotted for samples in each century. Each line is plotted with two probability intervals 0.65 (shorter thicker line) and 0.99 (longer thinner line). The probability density functions are also added to the top of the bar.
8.3.2 Trends by phase
The relationship between animal NISP counts and chronological phases is similar to the previous one, but with a different indexing variable for the intercept. Instead of being indexed by century, the intercept \(\alpha\) is now indexed by phase (\({[ChrID]}\)). While a phase can span several centuries, this modification allows the model to estimate separate intercepts for each phase, providing more specific insight into potential differences in the outcome variable and the average probability of animals occurrences over time. The rest of the model structure and distributional assumptions remain the same as the previous model, which uses a betabinomial distribution to model overdispersion in the data and an intercept-only structure.
\[ A_{i} \sim BetaBinomial(NISP_{i}, \bar{p}_{i} , \phi_{i}) \]
\[ logit(\bar{p}_{i}) = \alpha_{[ChrID]} \]
\[ \alpha_{[ChrID]} \sim Normal(0,1.5) \]
\[ \phi_{[ChrID]} \sim Exponential(1)+2 \]
As anticipated, the phase-level trends exhibit posterior distributions similar to the more detailed century-level results. However, the phase-level model offers an advantage - the credible intervals are considerably narrower than those of the previous model. Since each phase spans several centuries, except for the medieval phase that only includes the 11th century, there are more data points available to provide more reliable probabilities. Pigs’ NISP are at their highest point during the Roman phase, whereas they decrease in the Late Roman and Early Medieval phase, only to be slightly below the millennium average (of 0.37) in the 11th century. The predicted millennium average is comparable to the previous century-level model (0.36), with a small difference which can be due to the grouping of the observations. The diachronical means for the other animals are also very comparable to the other with a maximum difference of ±1 for each animal. Cattle’s NISP probabilities of occurrence increase throughout the millennium, although in the 11th century (marked as Ma on the plot) the 99% credible interval is larger, spanning from 0.15 to 0.26. Conversely, the century-level model for this period indicated a small increase towards the millennium average. As for the century-level model, the caprine are the second most represented animal in the millennium. Their probability of occurrence shows HDIs just below the millennium mean (of 0.27) in the Roman and Late Roman age, and positive HDIs values for the Early Medieval and Medieval phase, a situation not so different than the century-level model. The phase-level trends for domestic fowl and wild mammals show perhaps clearer trends than the century-level models. The HDI for domestic fowl, which is below the millennium mean of 0.05 in the Roman age, slightly increases in the later phases, with a peak in the 11th century. This peak is however less reliable than other phases, as the credible interval range is much larger. Wild mammals posterior probabilities on the other hand are stable around the mean (0.04) with a small increase in probability during the Early Medieval and Medieval phases. In addition to improved credible intervals as opposed to the century-level model, the probability intervals for the \(\phi\) shape parameter are also more reliable, as the PDFs are more narrow and most of the values are around the mean. Unfortunately, a smaller amount of samples for the 11th century produced more dispersed curves.
8.3.2.1 Community plot
After estimating the intercepts for each animal separately, they were compiled together to get a cohesive view of the probability of occurrence of the four animals. This approach is not ideal, but it can still be informative with regards to the patterns of animal exploitation in the archaeological record. The joint community plot allows to visually compare the probability of occurrence of each animal species over time and identify trends in their occurrence. For instance, the plot below shows how the predicted probability values for the three main domesticates—pigs, cattle and caprine vary across the four phases. Pigs have higher probabilities during the Roman phase, while already in the Late Roman period they get closer to other domesticates. It is in this period that the inter-sample variability of \(\phi\) seem to decrease, although the values are not very high. The precision \(\phi\) is instead very variable for edible wild mammals and domestic fowl, as it can be expected given that it is much rarer to find these animals in a sample for various reasons already discussed in the methods. The high variability in the precision of the medieval (11th century) intercept and \(\phi\) estimates can be explained by the lower number of samples. Overall, the most noticeable trends in the plot are the increased number of cattle in the early medieval period and the decreased number of pigs from the Late Roman phase. Modelling each animal species separately does not fully capture the complex interactions that exist among them, and future work should explore more advanced modelling approaches that account for the co-occurrence of multiple animal species in the same samples.
8.4 Context type
8.4.1 Pigs
To estimate the animals’ occurrences probability in each chronology and context, I used a betabinomial distribution to model overdispersion in the data. The \(A\) on the left side of the formula is the outcome variable—the animal NISP counts for each observation \(i\). This is an intercept-only model, where the intercept \(\alpha\) carries an interaction index \({[TCid]}\) as the model will provide estimates for each context type and chronology under examination. The \(\phi\) parameter indicates the precision in the Beta distribution, modelled by chronology and context type.
\[ A_{i} \sim BetaBinomial(NISP_{i}, \bar{p}_{i} , \phi_{i}) \]
\[ logit(\bar{p}_{i}) = \alpha_{[TCid]} \]
\[ \alpha_{[TCid]} \sim Normal(0,1.5) \]
\[ \phi_{[TCid]} \sim Exponential(1) + 2 \]
8.4.2 Cattle
8.4.3 Caprine
Warning: The `size` argument of `element_line()` is deprecated as of ggplot2 3.4.0.
ℹ Please use the `linewidth` argument instead.
8.4.4 Edible W. Mammals
8.4.5 Community plot
8.5 Macroregion
To estimate the animals’ occurrences probability in each chronology and macroregion, I used a betabinomial distribution to model overdispersion in the data. The \(A\) on the left side of the formula is the outcome variable—the animal NISP counts for each observation \(i\). This is an intercept-only model, where the intercept \(\alpha\) carries an interaction index \({[REGid]}\) as the model will provide estimates for each macroregion and chronology under examination. The \(\phi\) parameter indicates the precision in the Beta distribution, modelled by chronology and macroregion.
\[ A_{i} \sim BetaBinomial(NISP_{i}, \bar{p}_{i} , \phi_{i}) \]
\[ logit(\bar{p}_{i}) = \alpha_{[REGid]} \]
\[ \alpha_{[REGid]} \sim Normal(0,1.5) \]
\[ \phi_{[REGid]} \sim Exponential(1)+2 \]
8.5.1 Pigs
8.5.2 Cattle
8.5.3 Caprine
8.5.4 Edible W. Mammals
8.6 Geography
Animals distributions can vary across different geographical features. This research has considered plain, coast, hill and mountains as the most common geographical features in the Italian peninsula. Archaeological excavations where zooarchaeological remains have been analysed are located at low altitudes. Although as expected there are more mountain sites in Northern Italy, sampled sites are evenly placed on plains, coastlands and hills across the three Italian macroregions.
To estimate the animals’ occurrences probability in each chronology and geography, I used a betabinomial distribution to model overdispersion in the data. The \(A\) on the left side of the formula is the outcome variable—the animal NISP counts for each observation \(i\). This is an intercept-only model, where the intercept \(\alpha\) carries an interaction index \({[GEOid]}\) as the model will provide estimates for each geography and chronology under examination. The \(\phi\) parameter indicates the precision in the Beta distribution, modelled by chronology and geography.
\[ A_{i} \sim BetaBinomial(NISP_{i}, \bar{p}_{i} , \phi_{i}) \]
\[ logit(\bar{p}_{i}) = \alpha_{[GEOid]} \]
\[ \alpha_{[GEOid]} \sim Normal(0,1.5) \]
\[ \phi_{[GEOid]} \sim Exponential(1)+2 \]
8.6.1 Pigs
8.6.2 Cattle
8.6.3 Caprine
8.6.4 Edible W. Mammals
8.7 Altitude
The probability of occurrence of the most common faunal remains can be modelled against the elevation of sites in the four time periods under consideration. It is worth noting that the sites where the zooarchaeological remains have been found are not evenly distributed. In the Roman age, most sites investigated are located between 0 and 100 MSL, whereas after there is an increasing number of remains from sites between 100 and 400 MSL. Whether this reflects a real shift in settlement patterns is outside the aims of this study, but it might still be informative to visualise the different distribution of sites across elevations.
The proposed model to estimate the probability of occurrence as related to the altitude (the slope \(\beta\)) and chronology (\({[ChrID]}\)) uses a betabinomial distribution to model overdispersion in the data. The \(A\) on the left side of the formula is the outcome variable—the animal NISP counts for each observation \(i\). This is a simple intercept with slope model, where the intercept \(\alpha\) carries an index \({[ChrID]}\) as the model provides estimates for each chronology under examination. A single \(\phi\) parameter indicates the precision in the Beta distribution.
\[ A_{i} \sim BetaBinomial(NISP_{i}, \bar{p}_{i} , \phi_{i}) \]
\[ logit(\bar{p}_{i}) = \alpha_{[ChrID]} + \beta_{[ChrID]}\cdot Alt_{i} \]
\[ \alpha_{ChrID} \sim Normal(0,1.5) \]
\[ \beta_{ChrID} \sim Normal(0,1.5) \]
\[ \phi \sim Exponential(1)+2 \]
8.7.1 Pigs
8.7.2 Cattle
8.7.3 Caprine
8.7.4 Edible W. Mammals
8.7.5 Community plot
